{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# DRAG 脉冲\n",
    "*版权所有 (c) 2021 百度量子计算研究所，保留所有权利。*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 内容概要\n",
    "本教程将演示如何使用量脉实现一个采用 DRAG（Derivative Reduction by Adiabatic Gate）方案生成脉冲实现的 X 门，并且与传统的 $\\pi$ 脉冲方案进行对比，本教程的大纲如下：\n",
    "\n",
    "- 背景介绍\n",
    "- 准备工作\n",
    "- 定义 DRAG 波形\n",
    "- 量脉实现\n",
    "- 总结"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 背景介绍\n",
    "\n",
    "由于超导电路不是理想的二能级系统，因此需要考虑能级泄漏引入的误差。对于弱失谐量子比特，能量泄漏到第三能级会使量子比特的状态脱离计算空间。为了克服这个问题，研究人员提出了 DRAG 方案 \\[1\\]，它通过修正驱动脉冲的波形来消除能级泄露引入的误差。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 准备工作\n",
    "成功安装量脉后，可以按照本教程运行下面的量脉程序。要运行此教程，您需要从量脉（Quanlse）和其它常用的 Python 库导入以下包："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import the Hamiltonian module\n",
    "from Quanlse.QHamiltonian import QHamiltonian as QHam\n",
    "\n",
    "# Import the function for calculating infidelity\n",
    "from Quanlse.Utils.Infidelity import unitaryInfidelity\n",
    "\n",
    "# Import related operators\n",
    "from Quanlse.QOperator import driveX, driveY, driveZ, number, duff\n",
    "\n",
    "# Import waveforms and functions used to process the waveforms' data\n",
    "from Quanlse.QWaveform import gaussian, dragY1, QJob, QJobList, QWaveform\n",
    "\n",
    "# Import simulator interface for Quanlse Cloud Service\n",
    "from Quanlse.remoteSimulator import remoteSimulatorRunHamiltonian as runHamiltonian\n",
    "\n",
    "# Import matplotlib for graphing purposes\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Import numpy and scipy\n",
    "import numpy as np\n",
    "import math\n",
    "from scipy import integrate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "为了使用量脉云服务，我们需要获得一个 token 接入云端。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Import Define class and set the token for cloud service\n",
    "# Please visit http://quantum-hub.baidu.com\n",
    "from Quanlse import Define\n",
    "Define.hubToken = \"\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 定义 DRAG 波形"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "我们考虑一个 transmon 量子比特，通常通过向 XY 通道施加微波脉冲实现高保真度的 X 门。将量子比特第一激发态的能量和第二激发态能量分别记为 $\\omega_1$ 和 $\\omega_2$，驱动频率记为 $\\omega_d$。利用旋转波近似，该系统的哈密顿量可以写成 \\[2\\]：\n",
    "\n",
    "$$ \n",
    "\\hat H_R / \\hbar = \\delta_1 |1\\rangle \\langle 1|+\\delta_2 |2\\rangle \\langle 2|+\\frac{\\alpha_q}{2}\\hat a^{\\dagger}\\hat a^{\\dagger}\\hat a \\hat a+\\frac{\\varepsilon_x(t)}{2}\n",
    "\\left[ \\hat{a}^\\dagger + \\hat{a} \\right]+\\frac{\\varepsilon_y(t)}{2}\n",
    "i \\left[\\hat{a}^\\dagger - \\hat{a}\\right]\n",
    ",\n",
    "$$\n",
    "\n",
    "其中 $\\alpha_q = \\omega_2 -2\\omega_1$ 是系统的失谐强度。$\\delta_1 = \\omega_1-\\omega_d$ 和 $\\delta_2 = \\omega_2-\\omega_d$ 是量子比特跃迁频率相对于驱动频率的失调。此外 $\\varepsilon_x(t)$ 和 $\\varepsilon_y(t)$ 是施加到 XY 通道的脉冲函数。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "在理想情况下，我们可以忽略 transmon 量子比特的更高能级。我们通常会施加与量子比特频率相等的脉冲，即 $\\delta _1$ 设为零。欲想得到该 X 通道的脉冲，我们直接求解方程：\n",
    "$$\n",
    "\\int_0^{t_g}\\varepsilon_x(t)dt=\\theta. \n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "对于高斯波形 $\\varepsilon_G=Ae^{(t-\\tau)^2 /2\\sigma^2}-B$，上述待求解方程变为 $\\int_0^{t_g}\\varepsilon_G(t)dt=\\theta$，于是得到绕 x 轴转动 $\\theta$ 所对应的脉冲振幅 $A$：\n",
    "$$\n",
    "A=\\theta/\\left( \\int_0^{t_g}e^{-(t-\\tau)^2/2\\sigma^2}dt-t_ge^{-\\tau^2/2\\sigma^2} \\right),\n",
    "$$\n",
    "\n",
    "$$\n",
    "B=Ae^{-\\tau^2/2\\sigma^2}.\n",
    "$$\n",
    "在上面的等式中，$A$ 表示实现旋转量子门所需的脉冲振幅；而 $B$ 使得在开始时刻和结束时刻脉冲的振幅为零。\n",
    "\n",
    "在下面的代码中，我们首先设置系统的在布洛赫球上的旋转角度以及失谐性。然后，我们定义计算高斯波形脉冲值的函数（量脉提供了常用[波形](https://quanlse.baidu.com/api/Quanlse/Quanlse.QWaveform.html)的函数）。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "theta_x = np.pi # The angle of rotation\n",
    "Delta = -0.4 * 2 * np.pi # The anharmonicity in GHz\n",
    "\n",
    "# Calculate the parameters\n",
    "def intTheta(tg):\n",
    "    y = integrate.quad(gaussian(tg, 1, tg / 2, tg / 4), 0, tg)\n",
    "    return y[0]\n",
    "\n",
    "def calAx(tg):\n",
    "    return theta_x / (intTheta(tg) - gaussian(tg, 1, tg / 2, tg / 4)(0) * tg)\n",
    "\n",
    "def calBx(tg):\n",
    "    return calAx(tg) * gaussian(tg, 1, tg / 2, tg / 4)(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "DRAG 脉冲的波形和失谐分别为：\n",
    "$$\n",
    "\\varepsilon_y(t) = -\\frac{\\dot {\\varepsilon_x}(t)}{\\alpha_q}, \n",
    "$$\n",
    "$$\n",
    "\\delta_1(t) = -\\frac{\\varepsilon_x^2(t)}{2\\alpha_q}.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "这里，我们根据上面的方程建立控制脉冲 $\\varepsilon_x(t)$ 和 $\\varepsilon_y(t)$，并设置驱动脉冲的失谐 $\\delta_1$："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define the control waveforms\n",
    "def epsilonX(t, params):\n",
    "    tg = params['tg']\n",
    "    a = calAx(tg)\n",
    "    b = calBx(tg)\n",
    "    return gaussian(tg, a, tg / 2, tg / 4)(t) - b\n",
    "    \n",
    "def epsilonY(t, params):\n",
    "    tg = params['tg']\n",
    "    a = calAx(tg)\n",
    "    return dragY1(tg, a, tg / 2, tg / 4)(t) / Delta\n",
    "\n",
    "# Set the drive detuning  \n",
    "def delta1(t, params):\n",
    "    tg = params['tg']\n",
    "    lamda = np.sqrt(2)\n",
    "    return - epsilonX(t, {\"tg\": tg}) ** 2 / 2 / Delta"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 量脉实现"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "量脉将模拟和优化所需的系统信息存储在哈密顿量的字典中。首先，我们通过实例化一个 `QHamiltonian` 的对象 `ham` 创建一个空的哈密顿量。为了做一个对比，我们创建了两个空的哈密顿量，`ham` 不使用 DRAG 方法，`hamDrag` 将会使用 DRAG 方法。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create empty Hamiltonians\n",
    "ham = QHam(subSysNum=1, sysLevel=3, dt=0.2222)\n",
    "hamDrag = QHam(subSysNum=1, sysLevel=3, dt=0.2222)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "对于这一特定任务，系统哈密顿量可分为四个部分：\n",
    "$$\n",
    "\\hat H_R = \\hat H_{\\rm drift} + \\hat H_{\\rm xctrl} + \\hat H_{\\rm yctrl}+ \\hat H_{\\rm freq} ,\n",
    "$$\n",
    "其中 $\\hat H_{\\rm drift}= \\alpha_q\\hat a^{\\dagger}\\hat a^{\\dagger}\\hat a \\hat a \\, / \\, 2$ 表示量子比特的失谐性，是超导量子比特的内在特性，且与时间无关，因而我们可以通过调用 `addDrift()` 来添加该项。算符 $\\hat{a}^{\\dagger}\\hat{a}^{\\dagger} \\hat{a} \\hat{a}$ 在量脉中用 `duff()` 来表示，它将系统的能级数量作为输入参数。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Add the anharmonic terms\n",
    "ham.addDrift(duff, onSubSys=0, coef=Delta / 2.0)\n",
    "hamDrag.addDrift(duff, onSubSys=0, coef=Delta / 2.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "下一步，通过调用 `addWave()` 添加控制项 $\\hat H_{\\rm xctrl}=(\\hat a +\\hat a^{\\dagger})/2$、$\\hat H_{\\rm yctrl}=i(\\hat a -\\hat a^{\\dagger})/2$ 和 $ \\hat H_{\\rm freq}=\\hat a^{\\dagger}\\hat a $。我们可以在 `QOperator` 中找到相应的算符。这里我们创建一个 `QJobList` 的对象 `JobList` 来统一添加控制波形。代码如下:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Gate times\n",
    "t = np.arange(2., 9., 0.5)\n",
    "\n",
    "# Intialize array index\n",
    "jobList = ham.createJobList()\n",
    "jobListDrag = hamDrag.createJobList()\n",
    "for tg in t:\n",
    "    jobWaves = jobList.createJob()\n",
    "    jobWavesDrag = jobListDrag.createJob()\n",
    "    # Add Gaussian Wave of X control on the qubit 0\n",
    "    paraArgs = {\"a\": -0.5 * 2.0 * np.pi}\n",
    "    # Add wave for the job list without DRAG pulses\n",
    "    jobWaves.addWave(driveX, 0, QWaveform(epsilonX, 0, tg, {\"tg\": tg}))\n",
    "    # Add wave for the job list with DRAG pulses\n",
    "    jobWavesDrag.addWave(driveX, 0, QWaveform(epsilonX, 0, tg, {\"tg\": tg}))\n",
    "    jobWavesDrag.addWave(driveY, 0, QWaveform(epsilonY, 0, tg, {\"tg\": tg}))\n",
    "    jobWavesDrag.addWave(driveZ, 0, QWaveform(delta1, 0, tg, {\"tg\": tg}))\n",
    "    # Append this job to the job list\n",
    "    jobList.addJob(jobWaves)\n",
    "    jobListDrag.addJob(jobWavesDrag)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "为了对 DRAG 脉冲方案和传统的 $\\pi$ 脉冲方案进行完整的对比，我们在不同的门时间下计算量子门的保真度。使用量脉可以非常高效地完成这一任务，我们提供的函数 `runHamiltonian()` 允许用户向云端提交批量任务。它返回一个包含详细结果的字典列表，酉矩阵存储在每个字典的键 `\"unitary\"` 下。\n",
    "\n",
    "在本地设备上进行模拟可能需要很长时间，量脉提供的云服务可以显著加快这一过程。要使用量脉云服务，用户可以在百度量易伏网站 [http://quantum-hub.baidu.com](http://quantum-hub.baidu.com) 上获取 token，然后使用 `remoteSimulator()` 模块中的函数将任务提交到量脉云服务上。\n",
    "\n",
    "下一步是计算我们刚刚获得的实际酉矩阵和理想 X 门之间的距离。量脉提供了 `unitaryInfidelity()` 的函数，该函数能够根据下式计算失真度：\n",
    "$$\n",
    "{\\rm infid} =1- \\frac{1}{2}\\left|{\\rm Tr}(\\hat{\\sigma}_x P(U))\\right|.\n",
    "$$\n",
    "其中投影演化 $P(U)$（$U$是系统的演化算符）是描述投影到由最低的两个能量本征态 $|0\\rangle$ 和 $|1\\rangle$ 构成的计算空间的演化；$\\hat{\\sigma}_x$ 是我们想要实现的目标门。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create the arrays for storing gate infidelities\n",
    "errorX = np.zeros(len(t))\n",
    "errorXDrag = np.zeros(len(t))\n",
    "\n",
    "# Submit the job lists to Quanlse Cloud Service\n",
    "result = runHamiltonian(ham, jobList=jobList)\n",
    "resultDrag = runHamiltonian(hamDrag, jobList=jobListDrag)\n",
    "errorX = []\n",
    "errorXDrag = []\n",
    "\n",
    "for index in range(len(t)):\n",
    "    errorX.append(unitaryInfidelity(np.array([[0, 1], [1, 0]], dtype=complex), result[index][\"unitary\"], 1))\n",
    "    errorXDrag.append(unitaryInfidelity(np.array([[0, 1], [1, 0]], dtype=complex), resultDrag[index][\"unitary\"], 1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "最后，我们可以使用 Matplotlib 库对结果进行分析和可视化。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.semilogy(t, errorXDrag, label='With DRAG', marker='.')\n",
    "plt.semilogy(t, errorX, label='Without DRAG', marker='.')\n",
    "\n",
    "plt.xlabel('Gate Time (ns)')\n",
    "plt.ylabel('Infidelity')\n",
    "plt.title('X Gate')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "如上所示，我们消除了大部分的能级泄露误差。蓝色（DRAG 优化）曲线显示 DRAG 脉冲显著地提升了 X 门的保真度。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 总结\n",
    "本教程使用量脉呈现了 DRAG 的优化方案与传统的 $\\pi$ 脉冲方案之间的对比。用户可以通过点击这个链接 [tutorial-drag.ipynb](https://github.com/baidu/Quanlse/blob/main/Tutorial/CN/tutorial-drag-cn.ipynb) 跳转到此 Jupyter Notebook 文档相应的 GitHub 页面并且运行这个程序。我们鼓励用户尝试不同于本教程的参数值以获得最佳结果。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 参考文献\n",
    "\n",
    "\\[1\\] [Motzoi, Felix, et al. \"Simple pulses for elimination of leakage in weakly nonlinear qubits.\" *Physical review letters* 103.11 (2009): 110501.](https://link.aps.org/doi/10.1103/PhysRevLett.103.110501)\n",
    "\n",
    "\\[2\\] [Krantz, Philip, et al. \"A quantum engineer's guide to superconducting qubits.\" *Applied Physics Reviews* 6.2 (2019): 021318.](https://aip.scitation.org/doi/abs/10.1063/1.5089550)"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "2a31ed61199c5c13a03065ecec963a63da8631d96d1b9e695dac4715cb4eadb9"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
